Se usar a origem o arquivo RDA, já considera o merge feito. Caso contrário ele vai considerar como origem os arquivos CSV e vai gerar o rda a partir deles:

origem_dados_arq_rda = TRUE

Primeiro vamos abrir os dados e fazer o merge com a base de indicadores AFD (Adequação da formação docente):

if (origem_dados_arq_rda) {
  # Os dados estão na variável df
  load('data/complete_data.rda')
} else {
  # Lê a base de dados e notas do saeb
  load('microdados_nota_saeb_pe.rda')
  
  # Aqui temos 2039 registros
  df = df_microdados_nota_saeb
  df$NU_ANO_CENSO = as.factor(df$NU_ANO_CENSO)
  
  # Lê as bases de indicadores AFD (adequação da formação docente) e faz o merge
  # A base AFD está completa, continuamos com 2039 registros
  afd_2017 = read.csv2('dados/AFD_ESCOLAS_2017.csv', na.strings=c('--'), dec=',')
  afd_2019 = read.csv2('dados/AFD_ESCOLAS_2019.csv', na.strings=c('--'), dec=',')
  afd_2021 = read.csv2('dados/AFD_ESCOLAS_2021.csv', na.strings=c('--'), dec=',')
  
  afd_2017$NU_ANO_CENSO = 2017
  afd_2019$NU_ANO_CENSO = 2019
  afd_2021$NU_ANO_CENSO = 2021
  
  afd = rbind(afd_2017, afd_2019, afd_2021)
  
  df = merge(df, afd, by=c('CO_ENTIDADE', 'NU_ANO_CENSO'))
  
  # Lê as bases de indicadores ATU (média de alunos por turma) e faz o merge
  # A base está completa, continuamos com 2039 registros
  atu_2017 = read.csv2('dados/ATU_ESCOLAS_2017.csv', na.strings=c('--'), dec=',')
  atu_2019 = read.csv2('dados/ATU_ESCOLAS_2019.csv', na.strings=c('--'), dec=',')
  atu_2021 = read.csv2('dados/ATU_ESCOLAS_2021.csv', na.strings=c('--'), dec=',')
  
  atu_2017$NU_ANO_CENSO = 2017
  atu_2019$NU_ANO_CENSO = 2019
  atu_2021$NU_ANO_CENSO = 2021
  
  atu = rbind(atu_2017, atu_2019, atu_2021)
  
  df = merge(df, atu, by=c('CO_ENTIDADE', 'NU_ANO_CENSO'))
  
  # Lê as bases de indicadores IRD (média de regularidade do docente) e faz o merge
  # A base não está completa, faltam 29 registros (10 de 2017, 13 de 2019 e 6 de 2021 - passamos para 2010 registros aqui)
  ird_2017 = read.csv2('dados/IRD_ESCOLAS_2017.csv', na.strings=c('--'), dec=',')
  ird_2019 = read.csv2('dados/IRD_ESCOLAS_2019.csv', na.strings=c('--'), dec=',')
  ird_2021 = read.csv2('dados/IRD_ESCOLAS_2021.csv', na.strings=c('--'), dec=',')
  
  ird_2017$NU_ANO_CENSO = 2017
  ird_2019$NU_ANO_CENSO = 2019
  ird_2021$NU_ANO_CENSO = 2021
  
  ird = rbind(ird_2017, ird_2019, ird_2021)
  
  df = merge(df, ird, by=c('CO_ENTIDADE', 'NU_ANO_CENSO'), all.x = TRUE)
  
  # Lê o indicador socioeconômico.
  # Aqui temos os dados a cada 4 anos. Vou pegar o de 2019 e considerar
  # como 2019 e 2021. Poderia também pegar o de 2015 e considerar como 2017,
  # mas a forma de cálculo mudou, então vou deixar apenas 2019/2021 mesemo
  
  inse_2019 = read.csv2('dados/INSE_2019_ESCOLAS.csv')
  inse_2021 = inse_2019
  
  inse_2019$NU_ANO_CENSO = 2019
  inse_2021$NU_ANO_CENSO = 2021
  
  inse = rbind(inse_2019, inse_2021)
  
  df = merge(df, inse, by=c('CO_ENTIDADE', 'NU_ANO_CENSO'), all.x = TRUE)
  
  save(df, file='microdados_e_indicadores.rda')
}

Começa filtrando os dados de 2017 e removendo os registros com 88888 (1 registro em cada ano)

df = df[df$NU_ANO_CENSO != 2017, ]
# Ver se ainda precisamos remover esses registros, pois na prática agora já não estamos mais usando QT_PROF_MONITORES
df = df[df$QT_PROF_MONITORES != 88888,]

Para comparar modelos sempre usando os mesmos registros, vamos remover os dados onde IRD = NULL:

(Obs: Se no final do processo não formos usar IRD, podemos comentar essa linha)

df = df[!is.na(df$IRD), ]

Colocando as variáveis que o Bruno encontrou para 17-19-21 nesse modelo e checando:

df$PROP_INTEGRAL = df$QT_MAT_MED_INT/df$QT_MAT_MED

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_MODULOS + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 1) + IN_LOCAL_FUNC_SALAS_OUTRA_ESC + IN_BANHEIRO_EI + IN_ALIMENTACAO + IN_EJA_MED + IN_REFEITORIO + IN_SECRETARIA + IN_COZINHA + IN_LIXO_SERVICO_COLETA + as.factor(TP_OCUPACAO_PREDIO_ESCOLAR) + as.factor(TP_REGULAMENTACAO), df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_MODULOS + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_LOCAL_FUNC_SALAS_OUTRA_ESC + IN_BANHEIRO_EI + IN_ALIMENTACAO + 
##     IN_EJA_MED + IN_REFEITORIO + IN_SECRETARIA + IN_COZINHA + 
##     IN_LIXO_SERVICO_COLETA + as.factor(TP_OCUPACAO_PREDIO_ESCOLAR) + 
##     as.factor(TP_REGULAMENTACAO), data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.24148 -0.30994 -0.02962  0.27817  2.96966 
## 
## Coefficients:
##                                         Estimate Std. Error t value Pr(>|t|)
## (Intercept)                             4.718656   0.534543   8.827  < 2e-16
## NU_ANO_CENSO2021                       -0.083742   0.025276  -3.313 0.000947
## QT_SALAS_UTILIZA_CLIMATIZADAS           0.024322   0.003201   7.599 5.55e-14
## IN_QUADRA_ESPORTES_COBERTA              0.126191   0.028185   4.477 8.20e-06
## IN_MODULOS                              0.009972   0.030654   0.325 0.744989
## IN_ORGAO_GREMIO_ESTUDANTIL              0.106014   0.026111   4.060 5.19e-05
## IN_ORGAO_OUTROS                        -0.157368   0.031655  -4.971 7.50e-07
## PROP_INTEGRAL                           0.519714   0.030945  16.795  < 2e-16
## QT_PROF_ALIMENTACAO                    -0.028462   0.007731  -3.682 0.000241
## log(QT_DESKTOP_ALUNO/100 + 1)           0.983956   0.139063   7.076 2.39e-12
## IN_LOCAL_FUNC_SALAS_OUTRA_ESC          -0.093903   0.047803  -1.964 0.049693
## IN_BANHEIRO_EI                          0.808269   0.157901   5.119 3.52e-07
## IN_ALIMENTACAO                         -0.438112   0.240640  -1.821 0.068887
## IN_EJA_MED                             -0.179057   0.032248  -5.553 3.39e-08
## IN_REFEITORIO                           0.045206   0.033826   1.336 0.181636
## IN_SECRETARIA                          -0.099951   0.031552  -3.168 0.001570
## IN_COZINHA                             -0.140383   0.057735  -2.431 0.015166
## IN_LIXO_SERVICO_COLETA                  0.409588   0.106947   3.830 0.000134
## as.factor(TP_OCUPACAO_PREDIO_ESCOLAR)2  0.012595   0.061710   0.204 0.838305
## as.factor(TP_OCUPACAO_PREDIO_ESCOLAR)3  0.158213   0.054556   2.900 0.003792
## as.factor(TP_REGULAMENTACAO)1           0.003314   0.465910   0.007 0.994325
## as.factor(TP_REGULAMENTACAO)2           0.368975   0.504390   0.732 0.464584
##                                           
## (Intercept)                            ***
## NU_ANO_CENSO2021                       ***
## QT_SALAS_UTILIZA_CLIMATIZADAS          ***
## IN_QUADRA_ESPORTES_COBERTA             ***
## IN_MODULOS                                
## IN_ORGAO_GREMIO_ESTUDANTIL             ***
## IN_ORGAO_OUTROS                        ***
## PROP_INTEGRAL                          ***
## QT_PROF_ALIMENTACAO                    ***
## log(QT_DESKTOP_ALUNO/100 + 1)          ***
## IN_LOCAL_FUNC_SALAS_OUTRA_ESC          *  
## IN_BANHEIRO_EI                         ***
## IN_ALIMENTACAO                         .  
## IN_EJA_MED                             ***
## IN_REFEITORIO                             
## IN_SECRETARIA                          ** 
## IN_COZINHA                             *  
## IN_LIXO_SERVICO_COLETA                 ***
## as.factor(TP_OCUPACAO_PREDIO_ESCOLAR)2    
## as.factor(TP_OCUPACAO_PREDIO_ESCOLAR)3 ** 
## as.factor(TP_REGULAMENTACAO)1             
## as.factor(TP_REGULAMENTACAO)2             
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.464 on 1347 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.4603, Adjusted R-squared:  0.4519 
## F-statistic: 54.71 on 21 and 1347 DF,  p-value: < 2.2e-16

Remove as que são diferentes de 0:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI +  IN_EJA_MED + IN_SECRETARIA + IN_COZINHA + IN_LIXO_SERVICO_COLETA, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + IN_COZINHA + 
##     IN_LIXO_SERVICO_COLETA, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.21382 -0.31023 -0.03029  0.27216  3.12770 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    4.367712   0.119741  36.476  < 2e-16 ***
## NU_ANO_CENSO2021              -0.084437   0.025300  -3.337 0.000868 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.025771   0.003183   8.096 1.25e-15 ***
## IN_QUADRA_ESPORTES_COBERTA     0.132478   0.028033   4.726 2.53e-06 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.109719   0.026154   4.195 2.90e-05 ***
## IN_ORGAO_OUTROS               -0.171121   0.031423  -5.446 6.12e-08 ***
## PROP_INTEGRAL                  0.506926   0.029528  17.168  < 2e-16 ***
## QT_PROF_ALIMENTACAO           -0.031582   0.007694  -4.105 4.29e-05 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.023388   0.135027   7.579 6.42e-14 ***
## IN_BANHEIRO_EI                 0.816718   0.158232   5.162 2.81e-07 ***
## IN_EJA_MED                    -0.180029   0.027434  -6.562 7.52e-11 ***
## IN_SECRETARIA                 -0.084685   0.029149  -2.905 0.003729 ** 
## IN_COZINHA                    -0.144168   0.057282  -2.517 0.011956 *  
## IN_LIXO_SERVICO_COLETA         0.359742   0.105739   3.402 0.000688 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4661 on 1356 degrees of freedom
## Multiple R-squared:  0.453,  Adjusted R-squared:  0.4478 
## F-statistic: 86.38 on 13 and 1356 DF,  p-value: < 2.2e-16

INDICADORES AFD

Agora vamos colocar os dados dos indicadores no modelo que já temos:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI +  IN_EJA_MED + IN_SECRETARIA + IN_COZINHA + IN_LIXO_SERVICO_COLETA + AFD_CAT_1 + AFD_CAT_2 + AFD_CAT_3 + AFD_CAT_4 + AFD_CAT_5, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + IN_COZINHA + 
##     IN_LIXO_SERVICO_COLETA + AFD_CAT_1 + AFD_CAT_2 + AFD_CAT_3 + 
##     AFD_CAT_4 + AFD_CAT_5, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.23011 -0.31527 -0.03138  0.28078  2.91626 
## 
## Coefficients: (1 not defined because of singularities)
##                                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    3.8194151  0.4374096   8.732  < 2e-16 ***
## NU_ANO_CENSO2021              -0.1199442  0.0260394  -4.606 4.49e-06 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.0249523  0.0031491   7.923 4.79e-15 ***
## IN_QUADRA_ESPORTES_COBERTA     0.1245751  0.0276998   4.497 7.47e-06 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.1091843  0.0258332   4.227 2.53e-05 ***
## IN_ORGAO_OUTROS               -0.1628263  0.0311406  -5.229 1.98e-07 ***
## PROP_INTEGRAL                  0.4667244  0.0298901  15.615  < 2e-16 ***
## QT_PROF_ALIMENTACAO           -0.0280687  0.0076366  -3.676 0.000247 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  0.9773189  0.1337383   7.308 4.63e-13 ***
## IN_BANHEIRO_EI                 0.7879741  0.1563317   5.040 5.27e-07 ***
## IN_EJA_MED                    -0.1834140  0.0271696  -6.751 2.18e-11 ***
## IN_SECRETARIA                 -0.0927158  0.0288277  -3.216 0.001330 ** 
## IN_COZINHA                    -0.1440451  0.0566003  -2.545 0.011040 *  
## IN_LIXO_SERVICO_COLETA         0.2863198  0.1054426   2.715 0.006704 ** 
## AFD_CAT_1                      0.0098710  0.0044082   2.239 0.025305 *  
## AFD_CAT_2                      0.0160033  0.0064559   2.479 0.013301 *  
## AFD_CAT_3                      0.0033846  0.0043404   0.780 0.435644    
## AFD_CAT_4                     -0.0004123  0.0047401  -0.087 0.930706    
## AFD_CAT_5                             NA         NA      NA       NA    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4601 on 1352 degrees of freedom
## Multiple R-squared:  0.4686, Adjusted R-squared:  0.4619 
## F-statistic: 70.13 on 17 and 1352 DF,  p-value: < 2.2e-16

Com esse resultado, vamos manter apenas os indicadores 1 e 2 do AFD:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI +  IN_EJA_MED + IN_SECRETARIA + IN_COZINHA + IN_LIXO_SERVICO_COLETA + AFD_CAT_1 + AFD_CAT_2, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + QT_PROF_ALIMENTACAO + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + IN_COZINHA + 
##     IN_LIXO_SERVICO_COLETA + AFD_CAT_1 + AFD_CAT_2, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.20964 -0.31394 -0.02707  0.27972  2.92483 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    4.105936   0.126786  32.385  < 2e-16 ***
## NU_ANO_CENSO2021              -0.109422   0.025525  -4.287 1.94e-05 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.024751   0.003149   7.860 7.81e-15 ***
## IN_QUADRA_ESPORTES_COBERTA     0.125278   0.027718   4.520 6.73e-06 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.109836   0.025835   4.251 2.27e-05 ***
## IN_ORGAO_OUTROS               -0.164125   0.031136  -5.271 1.58e-07 ***
## PROP_INTEGRAL                  0.470141   0.029841  15.755  < 2e-16 ***
## QT_PROF_ALIMENTACAO           -0.028074   0.007638  -3.676 0.000246 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  0.973038   0.133751   7.275 5.85e-13 ***
## IN_BANHEIRO_EI                 0.789575   0.156445   5.047 5.10e-07 ***
## IN_EJA_MED                    -0.180232   0.027135  -6.642 4.47e-11 ***
## IN_SECRETARIA                 -0.093111   0.028834  -3.229 0.001271 ** 
## IN_COZINHA                    -0.144581   0.056639  -2.553 0.010799 *  
## IN_LIXO_SERVICO_COLETA         0.279648   0.105312   2.655 0.008014 ** 
## AFD_CAT_1                      0.007147   0.001262   5.663 1.81e-08 ***
## AFD_CAT_2                      0.011890   0.004973   2.391 0.016955 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4604 on 1354 degrees of freedom
## Multiple R-squared:  0.467,  Adjusted R-squared:  0.4611 
## F-statistic: 79.09 on 15 and 1354 DF,  p-value: < 2.2e-16

Estamos com muitas variáveis no modelo. Vamos fazer um drop1 pra ver como fica o R2

# Essa função auxiliar remove Y de X. 
# Se X for um data.frame, Y é interpretado como as colunas
# e o retorno é o data frame X sem as colunas
# Caso contrário, retorna x-y
remove_y_de_x = function(x, y) {
  if (class(x) == 'data.frame') {
    for (nome_variavel in y) {
      df[, nome_variavel] = NULL
    }
    return (df)
  }
  return( setdiff(x, y) )
}

printa_r2_drop1 = function(df, variaveis_interesse) {
  df_analise = df[, variaveis_interesse]
  todas_variaveis_explicativas = setdiff(variaveis_interesse, c('NOTA_SAEB'))
  
  r2 = c()
  
  for (var in todas_variaveis_explicativas) {
    variaveis_explicativas_menos_uma = remove_y_de_x(todas_variaveis_explicativas, var)
    df_temp_modelo = df_analise[, c(variaveis_explicativas_menos_uma, 'NOTA_SAEB')]
    modelo = lm(NOTA_SAEB ~ ., df_temp_modelo)
    r2 = c(r2, 100*summary(modelo)$r.squared)
  }
  df_analise_r2 = data.frame(todas_variaveis_explicativas, r2)
  df_analise_r2 = df_analise_r2[order(-df_analise_r2$r2),]
  
  print(df_analise_r2)
}
variaveis_interesse = c('NOTA_SAEB', 'NU_ANO_CENSO', 'QT_SALAS_UTILIZA_CLIMATIZADAS', 'IN_QUADRA_ESPORTES_COBERTA', 'IN_ORGAO_GREMIO_ESTUDANTIL', 'IN_ORGAO_OUTROS', 'PROP_INTEGRAL', 'QT_PROF_ALIMENTACAO', 'QT_DESKTOP_ALUNO', 'IN_BANHEIRO_EI', 'IN_EJA_MED', 'IN_SECRETARIA', 'IN_COZINHA', 'IN_LIXO_SERVICO_COLETA', 'AFD_CAT_1', 'AFD_CAT_2')
printa_r2_drop1(df, variaveis_interesse)
##     todas_variaveis_explicativas       r2
## 15                     AFD_CAT_2 46.41417
## 12                    IN_COZINHA 46.39845
## 13        IN_LIXO_SERVICO_COLETA 46.35360
## 11                 IN_SECRETARIA 46.21756
## 7            QT_PROF_ALIMENTACAO 46.12696
## 1                   NU_ANO_CENSO 45.90762
## 4     IN_ORGAO_GREMIO_ESTUDANTIL 45.90350
## 3     IN_QUADRA_ESPORTES_COBERTA 45.82937
## 9                 IN_BANHEIRO_EI 45.61072
## 5                IN_ORGAO_OUTROS 45.60024
## 14                     AFD_CAT_1 45.38541
## 10                    IN_EJA_MED 44.85330
## 8               QT_DESKTOP_ALUNO 44.61844
## 2  QT_SALAS_UTILIZA_CLIMATIZADAS 44.04608
## 6                  PROP_INTEGRAL 36.64310

De acordo com esses dados, vou remover variáveis de forma que o R2 fique acima dos 46 (que é o valor que estava antes de colocar indicador):

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI +  IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1, 
##     data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.20594 -0.31127 -0.03494  0.27686  3.03728 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    4.208073   0.070593  59.611  < 2e-16 ***
## NU_ANO_CENSO2021              -0.116104   0.025781  -4.503 7.26e-06 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.024239   0.003172   7.641 4.05e-14 ***
## IN_QUADRA_ESPORTES_COBERTA     0.126073   0.028028   4.498 7.44e-06 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.108311   0.026090   4.151 3.51e-05 ***
## IN_ORGAO_OUTROS               -0.164946   0.031302  -5.270 1.59e-07 ***
## PROP_INTEGRAL                  0.468051   0.030049  15.576  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.044651   0.134603   7.761 1.65e-14 ***
## IN_BANHEIRO_EI                 0.752954   0.158064   4.764 2.11e-06 ***
## IN_EJA_MED                    -0.199923   0.026681  -7.493 1.21e-13 ***
## IN_SECRETARIA                 -0.134165   0.027740  -4.837 1.47e-06 ***
## AFD_CAT_1                      0.007564   0.001261   5.998 2.55e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4658 on 1358 degrees of freedom
## Multiple R-squared:  0.4529, Adjusted R-squared:  0.4485 
## F-statistic: 102.2 on 11 and 1358 DF,  p-value: < 2.2e-16
plot(modelo)

Nesse ponto temos 11 variáveis.

INDICADORES ATU

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + ATU_TOTAL, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + 
##     ATU_TOTAL, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.19303 -0.30713 -0.03324  0.27822  3.00894 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    4.081541   0.096914  42.115  < 2e-16 ***
## NU_ANO_CENSO2021              -0.111242   0.025883  -4.298 1.85e-05 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.024070   0.003171   7.592 5.84e-14 ***
## IN_QUADRA_ESPORTES_COBERTA     0.122807   0.028054   4.378 1.29e-05 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.104274   0.026151   3.987 7.04e-05 ***
## IN_ORGAO_OUTROS               -0.160706   0.031351  -5.126 3.39e-07 ***
## PROP_INTEGRAL                  0.450171   0.031455  14.311  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.038589   0.134511   7.721 2.22e-14 ***
## IN_BANHEIRO_EI                 0.772141   0.158233   4.880 1.19e-06 ***
## IN_EJA_MED                    -0.196968   0.026700  -7.377 2.81e-13 ***
## IN_SECRETARIA                 -0.131565   0.027746  -4.742 2.34e-06 ***
## AFD_CAT_1                      0.007070   0.001286   5.497 4.61e-08 ***
## ATU_TOTAL                      0.004549   0.002390   1.904   0.0572 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4653 on 1357 degrees of freedom
## Multiple R-squared:  0.4544, Adjusted R-squared:  0.4496 
## F-statistic: 94.18 on 12 and 1357 DF,  p-value: < 2.2e-16

Se colocarmos todos os coeficientes no modelo de ATU (Total, série 1, série 2 e série 3), nenhum será estatisticamente diferente de 0. Colocando apenas ATU_TOTAL (a média dos 3), ele é estatisticamente diferente de 0 pra p < 0.1 apenas. Na verdade, se desconsiderarmos o merge da base IRD que faz reduzir 19 registros de 2019 e 2021, o coeficiente até fica estatisticamente dif de 0 com p < 0.05. É muito frágil, não vale a pena considerar essa variável. Vamos voltar pro modelo que tínhamos antes:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1, 
##     data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.20594 -0.31127 -0.03494  0.27686  3.03728 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    4.208073   0.070593  59.611  < 2e-16 ***
## NU_ANO_CENSO2021              -0.116104   0.025781  -4.503 7.26e-06 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.024239   0.003172   7.641 4.05e-14 ***
## IN_QUADRA_ESPORTES_COBERTA     0.126073   0.028028   4.498 7.44e-06 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.108311   0.026090   4.151 3.51e-05 ***
## IN_ORGAO_OUTROS               -0.164946   0.031302  -5.270 1.59e-07 ***
## PROP_INTEGRAL                  0.468051   0.030049  15.576  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.044651   0.134603   7.761 1.65e-14 ***
## IN_BANHEIRO_EI                 0.752954   0.158064   4.764 2.11e-06 ***
## IN_EJA_MED                    -0.199923   0.026681  -7.493 1.21e-13 ***
## IN_SECRETARIA                 -0.134165   0.027740  -4.837 1.47e-06 ***
## AFD_CAT_1                      0.007564   0.001261   5.998 2.55e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4658 on 1358 degrees of freedom
## Multiple R-squared:  0.4529, Adjusted R-squared:  0.4485 
## F-statistic: 102.2 on 11 and 1358 DF,  p-value: < 2.2e-16

INDICADOR IRD

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + IRD, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + 
##     IRD, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.16755 -0.30998 -0.03231  0.27155  2.91655 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    3.798546   0.113015  33.611  < 2e-16 ***
## NU_ANO_CENSO2021              -0.121651   0.025619  -4.749 2.27e-06 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.022963   0.003161   7.265 6.29e-13 ***
## IN_QUADRA_ESPORTES_COBERTA     0.119351   0.027859   4.284 1.96e-05 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.101353   0.025941   3.907 9.80e-05 ***
## IN_ORGAO_OUTROS               -0.167114   0.031073  -5.378 8.86e-08 ***
## PROP_INTEGRAL                  0.445950   0.030208  14.763  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.045714   0.133607   7.827 1.00e-14 ***
## IN_BANHEIRO_EI                 0.731105   0.156965   4.658 3.51e-06 ***
## IN_EJA_MED                    -0.211900   0.026610  -7.963 3.52e-15 ***
## IN_SECRETARIA                 -0.124563   0.027613  -4.511 7.01e-06 ***
## AFD_CAT_1                      0.007711   0.001252   6.159 9.63e-10 ***
## IRD                            0.130772   0.028315   4.619 4.23e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4623 on 1357 degrees of freedom
## Multiple R-squared:  0.4614, Adjusted R-squared:  0.4566 
## F-statistic: 96.87 on 12 and 1357 DF,  p-value: < 2.2e-16
plot(modelo)

Uma observação aqui é que até então eu estava simulando os dados apenas com IFD, aí depois coloquei os dados de ATU e, agora, IRD. Quando coloquei IRD, a base mudou, perdendo 29 registros.

O coeficiente IRD é estatisticamente diferente de 0 e igual a 0.13. A diferença entre o primeiro e o terceiro quartil são de 0.6 pontos, o que significa uma variação na nota do saeb de 0.078, que é bem pouca. Mas se considerarmos o mínimo e o máximo de amplitude nós teríamos uma variação maior, de 0.35. Essa variável ajuda a diminuir a diferença de algun pontos extremos.

Nesse ponto vamos rodar mais um drop1:

variaveis_interesse = c('NOTA_SAEB', 'NU_ANO_CENSO', 'QT_SALAS_UTILIZA_CLIMATIZADAS', 'IN_QUADRA_ESPORTES_COBERTA', 'IN_ORGAO_GREMIO_ESTUDANTIL', 'IN_ORGAO_OUTROS', 'PROP_INTEGRAL', 'QT_DESKTOP_ALUNO', 'IN_BANHEIRO_EI', 'IN_EJA_MED', 'IN_SECRETARIA', 'AFD_CAT_1', 'IRD')

printa_r2_drop1(df, variaveis_interesse)
##     todas_variaveis_explicativas       r2
## 4     IN_ORGAO_GREMIO_ESTUDANTIL 45.51989
## 3     IN_QUADRA_ESPORTES_COBERTA 45.41934
## 10                 IN_SECRETARIA 45.34219
## 12                           IRD 45.27867
## 8                 IN_BANHEIRO_EI 45.26022
## 1                   NU_ANO_CENSO 45.24685
## 5                IN_ORGAO_OUTROS 45.06420
## 11                     AFD_CAT_1 44.65829
## 2  QT_SALAS_UTILIZA_CLIMATIZADAS 43.90299
## 7               QT_DESKTOP_ALUNO 43.70799
## 9                     IN_EJA_MED 43.59509
## 6                  PROP_INTEGRAL 37.29252

Poderíamos substituir o IN_ORGAO_GREMIO_ESTUDANTIL pelo IRD e ainda ter um R2 melhor do que o modelo anterior. Mas por enquanto vamos manter todas essas variáveis.

Indicador INSE

Apesar desse indicador não poder ser colocado num dos 4 eixos, fiz as contas com ele também.

Primeiro é importante ver como o R2 aumenta consideravelmente, indo pra mais de 55%:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + IRD + INSE_VALOR_ABSOLUTO + INSE_QTD_ALUNOS, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + 
##     IRD + INSE_VALOR_ABSOLUTO + INSE_QTD_ALUNOS, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.08958 -0.27875 -0.01538  0.23954  2.10726 
## 
## Coefficients:
##                                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    1.3032798  0.1821502   7.155 1.37e-12 ***
## NU_ANO_CENSO2021              -0.0950413  0.0233604  -4.068 5.01e-05 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.0176053  0.0028917   6.088 1.49e-09 ***
## IN_QUADRA_ESPORTES_COBERTA     0.0887522  0.0255353   3.476 0.000526 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.1538491  0.0237455   6.479 1.29e-10 ***
## IN_ORGAO_OUTROS               -0.1082466  0.0285096  -3.797 0.000153 ***
## PROP_INTEGRAL                  0.4277628  0.0284934  15.013  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  0.9783283  0.1213074   8.065 1.60e-15 ***
## IN_BANHEIRO_EI                 0.6134511  0.1425092   4.305 1.79e-05 ***
## IN_EJA_MED                    -0.1253033  0.0248107  -5.050 5.01e-07 ***
## IN_SECRETARIA                 -0.1113892  0.0250714  -4.443 9.60e-06 ***
## AFD_CAT_1                      0.0025292  0.0012093   2.091 0.036672 *  
## IRD                            0.1442621  0.0258070   5.590 2.74e-08 ***
## INSE_VALOR_ABSOLUTO            0.6401812  0.0379947  16.849  < 2e-16 ***
## INSE_QTD_ALUNOS               -0.0012993  0.0001931  -6.728 2.52e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4193 on 1354 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.5579, Adjusted R-squared:  0.5533 
## F-statistic:   122 on 14 and 1354 DF,  p-value: < 2.2e-16
plot(modelo)

Agora vamos fazer aquele drop1 pra ter uma ideia do impacto no R2:

variaveis_interesse = c('NOTA_SAEB', 'NU_ANO_CENSO', 'QT_SALAS_UTILIZA_CLIMATIZADAS', 'IN_QUADRA_ESPORTES_COBERTA', 'IN_ORGAO_GREMIO_ESTUDANTIL', 'IN_ORGAO_OUTROS', 'PROP_INTEGRAL', 'QT_DESKTOP_ALUNO', 'IN_BANHEIRO_EI', 'IN_EJA_MED', 'IN_SECRETARIA', 'AFD_CAT_1', 'IRD', 'INSE_VALOR_ABSOLUTO', 'INSE_QTD_ALUNOS')

printa_r2_drop1(df, variaveis_interesse)
##     todas_variaveis_explicativas       r2
## 11                     AFD_CAT_1 55.56792
## 3     IN_QUADRA_ESPORTES_COBERTA 55.31242
## 5                IN_ORGAO_OUTROS 55.27315
## 1                   NU_ANO_CENSO 55.16299
## 8                 IN_BANHEIRO_EI 55.08305
## 10                 IN_SECRETARIA 55.06887
## 9                     IN_EJA_MED 54.84952
## 12                           IRD 54.67028
## 2  QT_SALAS_UTILIZA_CLIMATIZADAS 54.38795
## 4     IN_ORGAO_GREMIO_ESTUDANTIL 54.31199
## 14               INSE_QTD_ALUNOS 54.27052
## 7               QT_DESKTOP_ALUNO 53.66461
## 6                  PROP_INTEGRAL 48.13710
## 13           INSE_VALOR_ABSOLUTO 46.50446

Com isso, vamos deixar apenas as 4 variáveis mais importantes:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + INSE_VALOR_ABSOLUTO + INSE_QTD_ALUNOS, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + INSE_VALOR_ABSOLUTO + INSE_QTD_ALUNOS, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.27091 -0.29377 -0.02263  0.24012  2.16355 
## 
## Coefficients:
##                                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    1.2339404  0.1539896   8.013 2.38e-15 ***
## NU_ANO_CENSO2021              -0.0684351  0.0245169  -2.791  0.00532 ** 
## PROP_INTEGRAL                  0.5582632  0.0276136  20.217  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.2408912  0.1274435   9.737  < 2e-16 ***
## INSE_VALOR_ABSOLUTO            0.7460133  0.0368960  20.219  < 2e-16 ***
## INSE_QTD_ALUNOS               -0.0009704  0.0002048  -4.738 2.38e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4519 on 1363 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.4829, Adjusted R-squared:  0.481 
## F-statistic: 254.6 on 5 and 1363 DF,  p-value: < 2.2e-16
plot(modelo)

Note que mesmo com tão poucas variáveis o R2 ainda está maior do que o que tínhamos antes quando havíamos encerrado o modelo com o indicador IRD.

É ainda mais interessante ver que apenas o INSE e o ensino integral juntos explicam cerca de 44% da variação nos dados:

modelo = lm(NOTA_SAEB ~ NU_ANO_CENSO + PROP_INTEGRAL + INSE_VALOR_ABSOLUTO, df)
summary(modelo)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + PROP_INTEGRAL + INSE_VALOR_ABSOLUTO, 
##     data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.14080 -0.30823 -0.04716  0.25690  2.18697 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)          1.16775    0.15839   7.373 2.89e-13 ***
## NU_ANO_CENSO2021    -0.07106    0.02550  -2.787  0.00539 ** 
## PROP_INTEGRAL        0.64036    0.02740  23.369  < 2e-16 ***
## INSE_VALOR_ABSOLUTO  0.74041    0.03651  20.277  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.47 on 1365 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.4398, Adjusted R-squared:  0.4386 
## F-statistic: 357.2 on 3 and 1365 DF,  p-value: < 2.2e-16

Incluindo duas variáveis ao modelo completo (seção indicador IRD): quantidade de profissionais psicólogos, quantidade de profissionais coordenarores.

modelo_01 = lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + IRD, df)


modelo_02 <- lm(NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + IRD+
                  QT_PROF_PSICOLOGO + QT_PROF_COORDENADOR, df)

nobs(modelo_01)
## [1] 1370
nobs(modelo_02)
## [1] 1370
summary(modelo_01)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + 
##     IRD, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.16755 -0.30998 -0.03231  0.27155  2.91655 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    3.798546   0.113015  33.611  < 2e-16 ***
## NU_ANO_CENSO2021              -0.121651   0.025619  -4.749 2.27e-06 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.022963   0.003161   7.265 6.29e-13 ***
## IN_QUADRA_ESPORTES_COBERTA     0.119351   0.027859   4.284 1.96e-05 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.101353   0.025941   3.907 9.80e-05 ***
## IN_ORGAO_OUTROS               -0.167114   0.031073  -5.378 8.86e-08 ***
## PROP_INTEGRAL                  0.445950   0.030208  14.763  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.045714   0.133607   7.827 1.00e-14 ***
## IN_BANHEIRO_EI                 0.731105   0.156965   4.658 3.51e-06 ***
## IN_EJA_MED                    -0.211900   0.026610  -7.963 3.52e-15 ***
## IN_SECRETARIA                 -0.124563   0.027613  -4.511 7.01e-06 ***
## AFD_CAT_1                      0.007711   0.001252   6.159 9.63e-10 ***
## IRD                            0.130772   0.028315   4.619 4.23e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4623 on 1357 degrees of freedom
## Multiple R-squared:  0.4614, Adjusted R-squared:  0.4566 
## F-statistic: 96.87 on 12 and 1357 DF,  p-value: < 2.2e-16
summary(modelo_02)
## 
## Call:
## lm(formula = NOTA_SAEB ~ NU_ANO_CENSO + QT_SALAS_UTILIZA_CLIMATIZADAS + 
##     IN_QUADRA_ESPORTES_COBERTA + IN_ORGAO_GREMIO_ESTUDANTIL + 
##     IN_ORGAO_OUTROS + PROP_INTEGRAL + log(QT_DESKTOP_ALUNO/100 + 
##     1) + IN_BANHEIRO_EI + IN_EJA_MED + IN_SECRETARIA + AFD_CAT_1 + 
##     IRD + QT_PROF_PSICOLOGO + QT_PROF_COORDENADOR, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.17950 -0.30857 -0.03287  0.26802  2.94413 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                    3.809654   0.112575  33.841  < 2e-16 ***
## NU_ANO_CENSO2021              -0.112861   0.025580  -4.412 1.10e-05 ***
## QT_SALAS_UTILIZA_CLIMATIZADAS  0.019031   0.003292   5.782 9.18e-09 ***
## IN_QUADRA_ESPORTES_COBERTA     0.116645   0.027720   4.208 2.75e-05 ***
## IN_ORGAO_GREMIO_ESTUDANTIL     0.100530   0.025804   3.896 0.000103 ***
## IN_ORGAO_OUTROS               -0.153498   0.031126  -4.931 9.17e-07 ***
## PROP_INTEGRAL                  0.464279   0.030454  15.245  < 2e-16 ***
## log(QT_DESKTOP_ALUNO/100 + 1)  1.037392   0.133038   7.798 1.25e-14 ***
## IN_BANHEIRO_EI                 0.551438   0.166146   3.319 0.000927 ***
## IN_EJA_MED                    -0.209581   0.026480  -7.915 5.12e-15 ***
## IN_SECRETARIA                 -0.121170   0.027493  -4.407 1.13e-05 ***
## AFD_CAT_1                      0.006908   0.001263   5.468 5.41e-08 ***
## IRD                            0.131971   0.028184   4.682 3.12e-06 ***
## QT_PROF_PSICOLOGO              0.190159   0.075389   2.522 0.011770 *  
## QT_PROF_COORDENADOR            0.038471   0.014232   2.703 0.006955 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4599 on 1355 degrees of freedom
## Multiple R-squared:  0.4679, Adjusted R-squared:  0.4624 
## F-statistic:  85.1 on 14 and 1355 DF,  p-value: < 2.2e-16
plot(modelo_01)

plot(modelo_02)

As duas variáveis parecem explicar parte do que a variável banheiro EI explicava. O R quadrado melhora minimamente, mas no gráfico da distância de Cook, percebe-se que as duas variáveis introduzem outliers influentes no modelo, sobretudo a variável relativa à quantidade de coordenadores. Não há perda de observações ao introduzir essas variáveis. A variável sobre quantidade de psicólogos pode explicar o Eixo de Inovação, que busca fazer a escola mais atrativa ao jovem, tentando entender o que eleva as taxas de desistência e o desinteresse pelo estudo.